Angular collapse of dipolar Bose-Einstein condensates 
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We explore the structure and dynamics of dipolar Bose-Einstein condensates (DBECs) near their 
threshold for instability. Near this threshold a DBEC may exhibit nontrivial, biconcave density 
distributions, which are associated with instability against collapse into "angular roton" modes. 
Here we discuss experimental signatures of these novel features. In the first, we infer local collapse 
of the DBEC from the experimental stability diagram. In the second, we explore the dynamics of 
collapse and find that a nontrivial angular distribution is a signature of the DBEC possessing a 
biconcave structure. 



I. INTRODUCTION 

The role of dipole-dipole interactions in the struc- 
ture and dynamics of quantum many-body systems is 
at the forefront of both theoretical and experimental 
research. The dipole-dipole interaction is long-ranged 
and anisotropic, in contrast to the isotropic s-wave con- 
tact interaction, and leads to novel physics in ultracold 
spinor 0, fermonic 0-0] and bosonic (8l-fl2| many- 
body systems. 

For example, the anisotropy of the dipole-dipole in- 
teraction leads to unusual properties of a Bose-Einstein 
condensate (BEC). The condensate will experience "mag- 
netostriction," wherein its aspect ratio does not match 
that of the trap in which it is confined [HJ]. Addition- 
ally, a dipolar BEC (DBEC) is expected to exhibit novel 
density profiles, for instance a biconcave shape where the 
maximum density occurs at theperiphery of the conden- 
sate rather than at its center Q, or else more elaborate 
densities for non-cylindrically symmetric traps In 
contrast to magnetostriction, these density profiles have 
not thus far been directly imaged in experiments since 
shapes in the dentiy profiles likely wash out upon expan- 
sion of the condensate. 

Unusual density distributions are, however, intimately 
related to a more clearly observable property of the gas, 
namely, its macroscopic collapse. Because the dipole- 
dipole interaction is always attractive for dipoles in a 
head-to-tail orientation, a DBEC will always become un- 
stable for a sufficiently large dipole density and collapse 
inward, just as is the case for atomic BECs with attrac- 
tive contact interactions. The novelty is that, since the 
dipole-dipole interaction is anisotropic, the density at 
which the instability occurs depends on the aspect ra- 
tio of the harmonic trapping potential that confines the 
dipoles. 

This instability has been probed experimentally in a 
DBEC of atomic 52 Cr [15j. In this experiment the con- 
densate was stabilized by tuning the s-wave scattering 
length to a sufficiently large value. Instability was then 
triggered by reducing this scattering length below a crit- 
ical value. Strikingly, the observed anisotropic density 
distribtion of the collapsed cloud was well-reproduced 
within the mean- field theory (T^ |. 

On the theory side, the detailed mechanism of the col- 



lapse has been discussed in [13, LL8|- The picture that 
emerges is that for sufficiently oblate traps the collapse 
of the DBEC occurs via local density fluctuations, rather 
than a global collapse to the trap's center. In this case the 
instability is driven by fluctuations into a soft "roton" - 
like mode, which is unique to condensates with dipolar 
interactions 8]. Depending on the trap's aspect ratio, the 
roton can have nodal surfaces in either the radial or angu- 
lar coordinate, for a cylindrically symmetric condensate, 
and these nodal patterns should dictate the details of the 
condensate's collapse. In particular, for condensates with 
a biconcave density profile, the roton decay mode should 
always have an angular nodal pattern. Measuring an an- 
gular collapse would then be an experimental signature 
of the biconcavity. 

Thus far there is no direct experimental evidence for 
the roton or for the local collapse in a DBEC, although we 
have argued that the data in Ref. [l5| supports the idea 
of a local collapse [18]. Further, Ref. [TtJ has explored 
circumstances of local collapse in various trap geometries, 
contrasting approches where the collapse is initiated by 
either rapid or else adiabatic changes in the scattering 
length. 

In this article we therefore tackle head-on the prospects 
for observing local collapse. After some preliminary re- 
marks in Sec. [Til we proceed in Sec. IIIII to analyze and 
extend the experimental result in Ref. [Hj]. We argue 
that for oblate traps the scattering length required to 
stabilize the condensate can be explained within mean- 
field theory, but only if the collapse occurs locally. Fur- 
ther, the distinction between local and global collapse 
becomes more clear if the number of dipoles is increased. 
In Sec. IIVI we develop an understanding of a more di- 
rect measurement of condensate collapse, following the 
experimental procedure of Ref. [l6j . which includes an 
expansion that allows for imaging of the cloud. We show 
that angular structure in the expanded image is a direct 
signature of biconcave structure. 



II. PRELIMINARIES 

Consider a dilute gas of dipolar entities each with 
dipole moment d polarized in the axial (z) direction. The 
interaction between two such entities interacting by both 
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the dipole-dipole and contact interactions is given, in cgs 
units, by 

,., ,s , 2 l-3cos 2 6> Anh 2 a s .. 
V(r-r') = rf 2 , + — — <S(r-r') (1) 
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where |r — r'| is the distance between the particles, 9 is 
the angle between r — r' and the z-axis, M is the particle 
mass and a s is the s-wave scattering length of the parti- 
cles. The contact interaction (second term in Eq. ([1])) is 
either repulsive (a s > 0) or attractive (a s < 0), regardless 
of the orientation of the particles. The dipole-dipole in- 
teraction (first term in Eq. (|T])), however, changes sign 
depending on the particle's orientation. Two dipoles 
aligned in the direction of their polarization (9 = 0) at- 
tract each other while two dipoles aligned orthogonal to 
this direction (9 — tt/2) repel each other. 

We consider such a gas confined by a harmonic po- 
tential of the form U(r) = ^Muj 2 (p 2 + X 2 z 2 ) where 
A = uj z /ujp is the trap aspect ratio, describing to what 
degree the trap is prolate (A < 1) or oblate (A > 1). The 
trapping potential introduces a zero-point contribution 
to the condensate energy which serves to stabilize the 
system. A gas without dipoles but possessing a small neg- 
ative scattering length proves stable for sufficiently low 
density at any trap aspect ratio. The negative scattering 
length at which the condensate goes unstable scales only 
weakly with trap aspect ratio. When the stability thresh- 
old is crossed, e.g. when the scattering length becomes 
sufficiently negative to destabilize the BEC, the conden- 
sate undergoes macroscopic collapse. For purely contact 
interactions, the mean-field potential of the condensate 
is directly proportional to the density of the condensate, 
so collapse occurs where the particle density is greatest, 
at the center of the trap fl9| . 

By contrast, the trap aspect ratio plays a decisive role 
in determining the stability of a DBEC. In a prolate trap, 
a DBEC behaves much like a BEC with attractive contact 
interactions. This geometry favors attraction between 
dipoles and will induce a global collapse to the center for 
a critical dipole-dipole interaction strength. Vice versa, 
a DBEC in an oblate trap might be expected to behave 
much like a BEC with repulsive contact interactions since 
the dipolar entities are predominately repulsive in this 
geometry. However, as shown in [9J, there exists a finite 
critical dipole-dipole interaction strength, for any aspect 
ratio, at which a DBEC becomes unstable. The mecha- 
nism for collapse in this large A regime, however, is very 
different than that of a DBEC in a prolate trap or of a 
BEC with purely contact interactions. 

In an oblate trap, the axial trapping frequency is large, 
which acts to suppress elongation in the trap center, ren- 
dering the global collapse unlikely. Instead, the dipoles in 
the condensate are expected to form local density max- 
ima whose spatial widths are on the order of the axial 
harmonic oscillator length a z — y/h/Muj z . Each such 
bunch of dipoles then elongates axially, leading to local 
collapse. These local density maxima are related to the 



softening of a roton mode, whose characteristic wave- 
length a z sets the scale of the local collapse [H, . 



III. LOCAL COLLAPSE: EVIDENCE FROM 
THE STABILITY DIAGRAM 

Thus far only one experiment has explored the stabil- 
ity of a DBEC as a function of the trap aspect ratio [ll| . 
Rather than tune the dipole moment to a critical value, 
the experiment instead artificially stabilized the conden- 
sate by introducing a positive s-wave scattering length 
via a Fano-Feshbach resonance. Upon reducing the scat- 
tering length below a critical value a cr it, the experiment 
was able to trigger collapse in the DBEC. 

The resulting experimental stability diagram (repro- 
duced from [I5|]) is presented in Figure U as a plot of 
the critical scattering length a cr ;t versus aspect ratio A. 
These results represent the measurement performed on 
a condensate of N = 2 x I0 4 52 Cr atoms. For prolate 
traps, a comparatively large scattering length is required 
to achieve stability. As A is increased, the zero-point 
energy in the axial direction stabilizes the DBEC, and 
stable condensates are possible with a smaller critical 
scattering length. 

This figure also shows the results of two alternative 
numerical calculations of critical scattering length. In 
one, the theoretical division between stable (shaded) 
and unstable (unshaded) regions of parameter space is 
determined by numerically solving the nonlocal Gross- 
Pitaevskii equation (CPE) using the potential in ([T]). A 
second approach, already employed as an approximation 
in the experimental paper, shows the dividing line be- 
tween the stable and unstable regions as a dashed line. 
This approximation posits a Gaussian ansatz wave func- 
tion (normalized to unity): 
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where a p and a z are the variational parameters and aho 



yJh/Mu), where u> — yuj^oj z is the geometric mean trap 
frequency. Using this ansatz, the Gross-Pitaevskii energy 
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\ip(r)\ 2 / ip*(r')V(r-T')ip(r')dr' 



dr, (3) 



where N is the condensate particle number, is calculated 
for a given a^o to determine whether the energy E[ip, if)*] 
has a minimum, and thus to determine if the conden- 
sate is energetically stable. The presence of a minimum, 
local or global, corresponds to the presence of a stable 
ground state. A key feature of the Gaussian trial wave 
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function is that it always places the maximum density in 
the condensate's center, i.e., it is incapable of describing 
local collapse. For prolate traps, the maximum density 
is in the center. In this case the Gaussian ansatz and the 
numerical solution to the GPE agree with each other on 
the critical scattering length, and they both are in good 
agreement with the experimental result. 

Care must be taken, however, using this approxima- 
tion for oblate condensates. This can be seen in the 
A > 1 region of the stability diagram in Figure [1] where 
the Gaussian ansatz predicts a lower critical scattering 
length than does the GPE. We attribute this difference 
to the ability of the GPE to model local collapse. In- 
deed, for larger aspect ratios we observe local collapse 
into roton-like modes, as we will discuss in the next sec- 
tion, and as have been reported in Ref. [ljj ■ Further, the 
experimental determination of a cr it tends to show better 
agreement with the GPE prediction than with that of the 
Gaussian ansatz. We interpret this as experimental sup- 
port for the occurence of local collapse, albeit somewhat 
indirect evidence. However, the roton modes involved in 
collapse might have either radial or angular nodal struc- 
ture. This experiment does not make this distinction. 

Within the uncertainty in the experiment, the data 
in Figure Q] discriminates between the two methods, but 
one may wish for a clearer discrimination. We there- 
fore consider cases where the atom number is increased. 
The critical scattering length a cr it is shown in Figure [5] 
for DBECs with atom numbers of N = 10 4 , 10 5 , and 
10 6 . For a given trap, increasing the number of dipoles 
increases the relative importance of the dipole-dipole in- 
teraction, which acts to further destabilize the conden- 
sate. Thus, as predicted by the GPE, a cr it increases with 
increasing atom number. Vice versa, the Gaussian ansatz 
predicts a more stable condensate with increasing atom 
number. The difference between the two theoretical ap- 
proaches could then be clearly distinguished in such an 
experiment. 

Although Figure [5] plots only the domain of aspect ra- 
tios < A < 10 3 , it is straightforward to obtain the 
stability thresholds in the A <C 1 and the A ^> 1 lim- 
its for fixed uj. These limits are usefully described in 
terms of the characteristic dipole length, given by (in 
CGS units) HI 
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In the limit of large dipole-dipole interactions, 
Nadd/aho 1, the interaction-dependent term in 

E[ijj,ijj*} dominates over the other terms since it scales 
with N. Thus, in this limit, the condensate is unstable 
if this interaction-dependent term is negative. Now, to 
treat the limits of very small and large A, we consider, re- 
spectively, the harmonic trapping potential in the limits 
lo z — > and uj p — > 0. 

Taking the limit lo z — > corresponds to an infinitely 
prolate, or cigar-shaped trap. In this geometry, the 
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FIG. 1: (color online). The stability diagram of a DBEC of 
N — 2 x 10 4 52 Cr atoms plotted as critical scattering length 
versus trap aspect ratio. The points show the experimental 
results of [15(|. the shaded regions show the results of solving 
the GPE exactly and the dashed line shows the results of the 
Gaussian ansatz. The theoretical methods disagree as trap 
aspect ratio A increases, and the exact results fit the experi- 
mental data with great accuracy. The pink (darker) regions 
are where biconcave structure is predicted on the condensate 
profile. 



dipolar mean-field term reduces to a simple coupling to 
the condensate density, d 2 J ij)*{pf) ^=^- ^{p')dv' = 

—4:Trh 2 add\^P(p)\ 2 /M because the dipole-dipole interac- 
tion reduces to a delta-function in p for this geome- 
try. Thus, in the quasi-lD geometry, the total mean- 
field term becomes Airh 2 (a s — add) IVK/ )! 2 - Similarly 
for the limit lo p — > 0, corresponding to an infinitely 
oblate trap, the dipolar mean-field term reduces to 
87rh 2 add\4>(z)\ 2 /M , giving a total mean-field term of 
4irh 2 (a s + 2add) \ip(z)\ 2 for this geometry. Thus, in the 
limit Nadd/aho ^> 1, we find that the DBEC is unsta- 
ble when a s < add for A <C I and the DBEC is unstable 
when a s < —2add for A > 1. These limits are indicated 
in Figure [21 and the mean-field calculations are detailed 
in Appendix [A"l 



IV. LOCAL COLLAPSE: EVIDENCE FROM 
THE COLLAPSED CLOUD 

To take a closer look at the nature of collapse, it is 
necessary to track the collapse itself as a function of 
time. This, too, has been achieved in the 52 Cr DBEC 
experiments, which are excellently reproduced by time- 
dependent mean-field theory [16j . These experiments 
have not, however, focused directly on observing conse- 
quences of local collapse. Here we discuss the prospects 
of making such a measurement. 

Briefly, in such an experiment the scattering length 
is altered from a value where the condensate is stable 
against collapse to a somewhat lower value a < a cr ; t . Af- 
ter this transition, the atoms begin their collapse into 
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FIG. 2: (color online). For a DBEC of atomic 52 Cr, this 
figure illustrates the critical s-wave scattering length (below 
which the DBEC is unstable) as a function of trap aspect ra- 
tio A in units of the Bohr radius for uj = 2tt x 700 Hz. The 
blue (dotted), black (dot-dashed) and red (solid) lines corre- 
spond to N = 10 4 , N = 10 5 and TV = 10 6 , respectively. The 
lines without crossed symbols are the results of the Gaussian 
ansatz and lines with crossed symbols are the results of the 
exact solution of the GPE. Notice how, as N is increased, the 
Gaussian ansatz predicts a more stable condensate while the 
exact solution predicts a less stable condensate. 



high density regions where three-body recombination 
takes over, ejecting atoms from the trap. The trap is 
generally released after some hold time, to expand the 
cloud for imaging. The resulting density patterns show 
intricate shapes and depend on details such as whether 
the passage from stable to unstable is adiabatic or dia- 
batic Fl7j ] . 



A. Modes of instability 

The underlying physics of the instability and collapse is 
determined by the softening of the roton modes. We com- 
pute these modes by solving the Bogoliubov de Gennes 
equations as in Ref. [loj . exploiting the cylindrical sym- 
metry of the system. Namely, we make the quasiparticle 
ansatz, 

V>(r, t) -> [i/)(p, z) + 5u{p, zJe^-^H- 

+8v*{p,z)e-^ mv -^A e-^ 4 (5) 

where uj is the quasiparticle energy, to is the projec- 
tion of the quasiparticle momentum onto the z-eods, \x 
is the chemical potential of the ground state ip(p, z) and 
5 <C 1 to ensure that the quasiparticles have small am- 
plitudes. By solving the BdG equations for various to 
quantum numbers, we determine whether the DBEC is 
dynamically stable or unstable by determining whether 



the quasiparticle energy is purely real or has a nonzero 
imaginary part, respectively [2fj| . 

Figure [3] illustrates the mode softening for a DBEC 
containing N = 10 4 52 Cr atoms at an aspect ratio A = 8. 
Plotted is the energy of the excitation as a function of the 
scattering length a s , labeled by its azimuthal angular mo- 
mentum quantum number to. The solid lines depict the 
real parts of these energies, while the symbols represent 
their imaginary parts. As a s diminishes, the energies of 
these modes drop to zero, and thereafter become purely 
imaginary. The first such transition, at a s ~ — 0.9ao, 
identifies the scattering length at which the DBEC is dy- 
namically unstable, since any small perturbation is capa- 
ble of exciting this mode, which then grows exponentially 
in time. Thus an unstable condensate quickly grows high- 
density peaks in regions defined by the antinodes of these 
modes. 

Figure [3] is a particular example illustrating the modes 
that contribute to instability at a particular aspect ratio 
A = 8. At this aspect ratio the condensate's density ex- 
hibits a biconcave shape, and so decay into angular rotons 
is expected. We stress that at all aspect ratios where the 
maximum density lies at the center rather than at the pe- 
riphery, the rotons responsible for instability are always 
to = modes that do not exhibit an angular structure. 
This connection is essential to connecting observed an- 
gular decay circumstantially to biconcave structure. 

Regardless of whether the roton is purely radial or an- 
gular in nature, it leads a DBEC to instability at a fixed 
length scale, as mentioned in Sec. [TTl having wavelength 
~ 2ira z . As the trap aspect ratio A is increased, the ra- 
tio of the axial to the radial harmonic oscillator lengths, 
a z /a pi is decreased, so more roton wavelengths can fit 
around the circumference of the condensate for larger A. 
For biconcave condensates, this results in angular rotons 
with larger to quantum number being responsible for in- 
stability for larger A, since more angular nodes can fit 
into the condensate in this regime. For N — 10 4 52 Cr 
atoms in a trap with A = 8 and uj — 2ir x 700 Hz, this 
mode happens to have to = 3. Indeed, the circumference 
of the region of maximum particle density in this bicon- 
cave condensate is ~ 67ra z , or three roton wavelengths. 



B. Numerics and the "ideal experiment" 

The mode that brings about the dynamical instability 
determines not only the scattering length at which the 
condensate will collapse, but also how the condensate will 
collapse as the stability threshold is crossed. Consider 
preparing a DBEC of N = 10 4 52 Cr atoms just above 
the stability threshold in a trap with aspect ratio A = 
8. These are the collective modes whose energies are 
shown in Figure [3J A small jump in scattering length to 
a value just below a cr j t would cause the condensate to go 
unstable by a macroscopic occupation of the to = 3 mode 
that has a nonzero imaginary energy at this scattering 
length. The density of the condensate during the collapse 
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FIG. 4: Collapsed condensate in a trap with mean frequency 
uj = 2n x 700 Hz and aspect ratio A = 8 after 10.5 ms. The 
perturbation for the collapse is controlled to have m = 3 
symmetry and the global phases as shown in the frames. Each 
frame corresponds to a different value of the initial phase a. 
The collapsed condensates are rotated by a/3, ensuring that 
the finite grid size does not influence the small length scale 
dynamics of the condensate collapse. 



FIG. 3: (color online). The real and imaginary parts of the 
low- lying BdG modes for a condensate of N = 10 4 52 Cr atoms 
with mean trap frequency oj = 2n x 700 Hz and trap aspect 
ratio A = 8, plotted as a function of the s-wave scattering 
length a s . The real parts are represented by solid lines and 
the imaginary parts, developing where the real parts go to 
zero, are represened by markers. The m — 3 mode, being the 
first to develop a nonzero imaginary energy, serves to define 
(Jcrit for this aspect ratio. 



would change, on a time scale r ~ 27r/Im[cj], as the atoms 
macroscopically occupy three clumps that self-attract in 
the z-direction. 

Decay of the condensate into a roton mode with m > 
requires breaking the condensate's initial cylindrical sym- 
metry by introducing fluctuations into the mode. In an 
experiment this is caused by thermal fluctuations, but in 
our calculation we must make this happen artificially. To 
do this, we seed the condensate wave function by adding 
to it a small contribution of the excited state mode: 

ip{r) -> </>(>, z) + 0.01e 2 ™ Q e 3 ^u 3 (p, z), (6) 

where e 3l< ^ describes the basic angular variation of the 
roton mode, and a is an additional phase that will de- 
termine the overall rotation of the collapsed condensate. 
In the absence of a seed like this, the numerical solu- 
tion remains at its unstable equilibrium for a time long 
compared to the natural lifetime 27r/Im[cj]. The appar- 
ent lifetime in this case is determined by the time be- 
fore roundoff error starts to affect the time evolution of 
the GPE. However, once the wave function is seeded as 
above, the decay occurs on the expected time scale. 

After the collapse is triggered, the condensate indeed 
forms the three clumps as expected, as seen in Figure [4] 
Shown is the density of particles, as viewed in the x-y 
plane, i.e., looking down from the axis of the dipoles' po- 
larization. Each peak was initially seeded by a density 
fluctuation at the antinode of the m = 3 excited state 
roton wave function. Thus the three peaks are uniformly 
equally spaced in angle, as befits the symmetery of the 
mode. An angular display of this sort would provide un- 



ambiguous evidence for nonlocal collapse. Moreover, the 
fact that the collapse occured in an angular coordinate 
provides indirect evidence for the biconcave structure of 
the initial state. 

While the relative positions of the three peaks in this 
experiment are well-defined by the symmetry of the roton 
mode, there is still an overall undetermined angle of rota- 
tion of the whole pattern. Numerically, this is set by the 
angle a in (J6)). Since the angular dependence of the con- 
densate density, with this wave function, is proportional 
to cos (3</5 + a), we expect that, if there is no unphysical 
dependence on the numerical grid, the collapse will occur 
rotated by an angle a/3 for any initial phase a. Indeed, 
we find that the collapse dynamics are unaffected by the 
grid, as is illustrated in Figure Here, we input the 
initial phases a = 7r/vl3 and a — 7r/v3 and find that 
the collapsed wave function is rotated by exactly these 
phases times 1/3. Although not shown here, simulations 
for other initial phases give the same results. Based on 
this ability to reproduce the same angular pattern, but 
rotated in a predictable way, we conclude that the under- 
lying Cartesian grid is adequate to describe this collapse. 



C. More realistic experiment 

Figure 2] illustrates the kind of clean angular distri- 
bution that might be expected in the ideal experiment, 
where an infinitesimal change in scattering length is pos- 
sible, and where only a single roton mode is excited. Thus 
far, neither of these circumstances is true in the 52 Cr 
experiments. Whereas Figure [3] shows the difference in 
scattering length at which modes with different m quan- 
tum numbers develop imaginary energies to be a fraction 
of a Bohr, experimental uncertainties in the Feshbach- 
tuncd scattering length of 52 Cr are ±2ao |l5(. Additon- 
ally, imaging of a 52 Cr DBEC was done after a time of 
free-expansion in this experiment, not in-trap as is de- 
scribed in the scenario above. We propose, with slight 
modification, an experiment similar to the one described 
by [lj| that presents us with the possibility of observ- 
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our condensate by adding excited modes (quasiparticles) 
with weights determined by the Bose-Einstein distribu- 
tion ElL 



(7) 



where rij is the number of particles occupying the quasi- 
particle state with energy uij , T is the temperature of the 
Bose gas, /j, is the chemical potential of the condensate 
and ks is the Boltzmann constant. 

Using the quasiparticle spectrum given by solving the 
BdG equations and a temperature of T = 100 nK, we 
then perturb our initial condensate by 



FIG. 5: (color online). Stability diagram for TV = 10 4 52 Cr 
atoms. The white region is dynamically unstable, while the 
darker regions are stable. The pink (darker) islands are where 
biconcave structure is found in the ground state of the con- 
densate. The mean trap frequency is Q = 2-k x 700 Hz for 
all aspect ratios A. The arrows illustrate the initial and final 
values of scattering length in the experiment proposed in the 
text. 



ing angular structure in the collapse and expansion of a 
DBEC. 

Instead of making a very small jump in scattering 
length across the stability threshold, consider making a 
jump of Aa s = -14a . For a 52 Cr DBEC with N = 10 4 
atoms, we numerically prepare, for Q ~ 2n x 700 Hz, 
a condensate in a trap with A = 2 and scattering length 
a s = 20ao and a condensate in a trap with A = 8 and scat- 
tering length a s = 10ao, where both scattering lengths 
are about 10ao above a cr it for their respective aspect ra- 
tios. We then ramp the scattering length from its initial 
value to its final value over a time period of 8 ms. These 
scattering length ramps are illustrated by the arrows in 
Figure [5] Although an 8 ms ramp time is not sufficiently 
slow to make the change completely adiabatic (the char- 
acteristic trap period is 2it/uj = 1.4 ms), it is sufficiently 
slow to allow a biconcave shape to form during the ramp. 
Once this ramp has been made, we hold the collapsing 
condensate in the trap for thoid = 2 ms and then turn 
off the trap to let the collapsed condensate propagate in 
free space. In an actual experiment, the expanded cloud 
could then be imaged to determine its density profile af- 
ter expansion. 

To ensure that we accurately simulate an experimen- 
tal scenerio and to break the cylindrical symmetry of 
the condensate in a physically consistent way, we seed 
the condensate prior to the time evolution with numer- 
ical noise. Although the non-projected GPE can only 
model the condensate dynamics at zero temperature, an 
experiment will unavoidably have a very small but fi- 
nite temperature present in the gas. This means that 
the condensate fraction will not be exactly unity, but 
something slightly less. To account for this, we correct 



j 

+<>,z)e-"^], (8) 

where {aij} are random numbers betwen and 1, rij is 
given by Eq. (J?} and u m j(p,z) and v^^p^z) are BdG 
modes with quantum number m and energy uij. Also, 
we include the factor \/l/N in the weighting because 
the condensate wavefunction ip(r) is normalized to unity 
instead of being normalized to N. We impose a cutoff on 
the sum in Eq. (@J of ojj < 2k B T, where T = 100 nK, 
simplifying the problem by ignoring higher energy modes 
that contribute little to the thermal excitations of the 
system. Indeed, setting T = 100 nK is an experimentally 
accessible temperature [22j ■ 

Additionally, because the condensate density becomes 
very large during the collapse process, a three-body loss 
term is required to accurately model the collapse and 
expansion dynamics [l5j . The rate constant for three- 
body recombination was experimentally determined to 
be L3 = 2 x 10~ 40 m 6 /s for 52 Cr. We account for this 
loss in our simulations by including the term —ihN(N — 
l)L 3 \tp(r)\ 4 /2 in the time-dependent GPE, given by 



+{N-1) J dr'y(r-r')|V'(r',t)| 2 
-A(iV-l)^|^r,t)| 4 }^(r,t), (9) 

where V(r — r') is given in Eq. (p}. 

Figure |6] illustrates the numerical time evolution of 
these condensates through the collapse and expansion de- 
scribed above. As before, these images represent density 
profiles in the x-y plane, as viewed from the polarization 
axis. The top four frames illustrate the collapse and ex- 
pansion of a condensate in a trap with aspect ratio A = 2, 
in which there is no biconcave shape and in which, con- 
sequently, there should be no collapse to angular roton 
modes. During its collapse, the condensate maintains its 
peak density in the center. After the trap is removed and 
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the gas is allowed to expand, its cylindrical symmetry is 
preserved. 

By contrast, the lower four panels of Figure [6] illustrate 
a representative time evolution for a condensate in a trap 
wtih aspect ratio A = 8. In this case, by 7.5 ms the con- 
densate has established its biconcave structure. When 
the condensate collapses, it does so into roton modes 
with angular nodal structure, leading to local collapse 
with angular nature. After the trap is turned off and 
the condensate expands for 4.5 ms, the angular struc- 
ture remains in the density of the expanded cloud. The 
collapse is clearly dominated by a roton with m = 3 in 
this simulation. However, because several angular modes 
are involved, the angular pattern no longer experiences 
pure m = 3 angular symmetry. Moreover, each mode ar- 
rives with a random initial phase, meaning that there is a 
random assymmetry due to the interference between the 
unstable modes. In the experiment this will imply non- 
repeatability of the observed density peaks from shot to 
shot. 

Nevertheless, once the angular pattern is established, 
its vestiges remain in the expanded cloud. In the final 
expanded picture, the clear break from cylindrical sym- 
metry indicates that the decay modes have angular de- 
pendence, hence that the condensate went through a bi- 
concave phase. 

We note that a collapse and expansion experiment was 
done on a DBEC of 52 Cr [23| . However, it did not probe 
the parameter regime for biconcave structure formation. 

The results of simulations that are very similar to the 
ones described above are presented in [lTj, where DBEC 
collapse is modeled in-trap and not through the expan- 
sion process and not including a three-body loss term in 
the simulation. Ref. pjj performs simulations of DBEC 
collapse for both adiabatic and non-adiabatic (instanta- 
neous) jumps in scattering length, and find very inter- 
esting results regarding the presence of global and local 
collapse in the condensate dynamics. For adiabatic col- 
lapse, where the change in scattering length is sufficiently 
slow to track the condensate across the roton softening in 
the BdG spectrum, they present a critical trap aspect ra- 
tio above which local collapse occurs. We confirm these 
results, but point out that while local collapse is very 
interesting (and can be evidence for the presence of the 
roton in these systems), its manifestation in a DBEC is 
much richer than has been discussed in previous work. A 
mapping of DBEC collapse via the experiment proposed 
above can determine not only whether collapse was global 
or local, but whether collapse was radial or angular and 
thus provide evidence for the underlying biconcave struc- 
ture. 

We point out that the experiment proposed above 
is just one of many experimental methods that would 
demonstrate the angular nature of DBEC collapse. Cer- 
tainly, taking data for a number of additional trap aspect 
ratios would assist in mapping out the regions where bi- 
concavity exists. Also, we expect that smaller and slower 
jumps in scattering length, which may be had with less 



uncertainty in the Feshbach-induced scattering length, 
would assist in understanding how and where collapse oc- 
curs. Slower ramping of the scattering length allows the 
condensate to be tracked more adiabatically and thus 
allows for collapse to begin when only one BdG mode 
has a nonzero imaginary energy, making the mapping of 
the collapse much more clear. Instantaneous or very fast 
jumps in scattering lengths across a biconcave region will 
miss this structure completely and thus result in a purely 
radial collapse. For an angular collapse to occur, the bi- 
concave structure must manifest itself in the condensate 
prior to collapse. 



V. CONCLUSION 

In conclusion, we have shown that in order to cor- 
rectly map out the stability of a DBEC, a computational 
method that is sensitive to the local nature of DBEC col- 
lapse must be used. Methods such as the gaussian ansatz 
that are not sensitive to such phenomenon will incorrectly 
predict the stability of the system. Also, we draw a con- 
nection between the BdG spectrum of a DBEC and the 
nature of the DBEC collapse. Not only can the BdG 
quasiparticles predict where a DBEC will collapse in pa- 
rameter space, they can also predict how a DBEC will 
collapse. For DBECs without biconcave structure, this 
collapse is purely radial while for DBECs with biconcave 
structure, this collapse has angular structure. Perform- 
ing collapse and expansion experiments on a 52 Cr DBEC 
can reveal this angular structure and thus provide an ex- 
perimental method for mapping biconcave structure in 
DBECs. 
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Appendix A: Calculation of the mean-field in 
reduced dimensions 

We consider the calculation of the mean-field poten- 
tial due to the dipole-dipole interaction in two different 
geometries, one with uj p — ¥ (a quasi-two dimensional 
(2D) geometry) and one with u> z — > (a quasi-one di- 
mensional (ID) geometry). In the quasi-2D geometry, we 
assume that the condensate wavefunction depends only 
on z and is homogeneous in the p-direction (or in the x- 
and y-directions) and in the quasi- ID geometry, we as- 
sume that the condensate wavefunction depends only on 
p (or on x and y) and is homogeneous in the z-direction. 
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x (cm) x-io x (cm) x-io # (cm) x-io % (cm) x-io 

FIG. 6: Collapse dynamics of DBECs, both in harmonic traps with mean frequency ui = 2n x 700 Hz, corresponding to the 
scattering length ramps illustrated in Figure [5] (a) A DBEC in a trap with aspect ratio A = 2, at t — ms the condensate has 
a 3 — 20ao, the scattering length is ramped down to a s = 6ao over 8 ms, the collapsed condensate is held in the trap for thoid = 2 
ms and then expanded until t — 16 ms. The collapse and expansion is purely radial, (b) A DBEC in a trap with aspect ratio 
A = 8, at t — ms the condensate has a s = 10ao, the scattering length is ramped down to a s — — 4ao over 8 ms, the collapsed 
condensate is held in the trap for thoid = 2 ms and then expanded until t = 14.5 ms. The condensate becomes biconcave during 
the ramp in scattering length and thus collapses with angular structure, preserving an angular character during expansion. 



We begin with the expression for the dipole-dipole in- 
teraction potential in momentum-space |24| 



V dd (k) 



-in 



d 2 (3 cos 2 k - 1) 



where 6*k is the angle between the direction of the dipole 
polarization (z or fc z , for the DBEC we are considering) 
and the vector k. Using this momentum-space repre- 
sentation, the coordinate-space mean-field potential due 
to the dipole-dipole interaction is given by the convolu- 
tion of Vdd(k) with the condensate density in momentum- 
space, n(k), 



UM^T- 1 Vrf d (k)n(k) 



(A2) 



where J 7-1 is the inverse Fourier Transform operator. 
First, consider the quasi-2D geometry, in which the con- 
densate density is homogeneous is x and y. The con- 
densate density in momentum-space is then given by the 
Fourier Transform, 

^2D(k) = T [t!2d(z)] = n 2 i)(k z )S(k x )S(ky). (A3) 

Substituting this result and Eq. (|A1[) into Eq. (|A2[) and 



writing cos 2 0k = k 2 / '(fc 2 + fc 2 + fc 2 ) gives an expression 
for the mean-field potential in the quasi-2D geometry, 



(Al) U 2 f(r)=T-i 



3 \ /c? ~\~ k„. ~\~ /j? 



1 ) n 2 B(k z )S(k x )S(k y ) 
(A4) 



The operation of the inverse Fourier Transform on this 
momentum-space function gives 

Uif(r) = 8 -^d 2 mz)\ 2 = ^^|^(z)| 2 , (A5) 

where |?/>(2:)| 2 is the coordinate-space condensate density 
in the quasi-2D geometry. We carry out the same calcu- 
lation for the quasi- ID geometry, where the condensate 
density in momentum-space is given by 



nm(k) = T [ni D (x, y)} = nm(k x ,k v )8(k z ) 
Substituting this function into Eq. (|A2|) gives 



(A6) 



4tt 



jnh add 
M 



mp)\ 2 , (A7) 



where IV'Cp)! 2 is the coordinate-space condensate density 
in the quasi- ID geometry, written in terms of p instead 
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of x and y. 
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